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Abstract 

A semiconductor-based terahertz-detector strategy, exploiting a bound-to-bound-to-continuum 
architecture, is presented and investigated. In particular, a ladder of equidistant energy levels is 
employed, whose step is tuned to the desired detection frequency and allows for sequential multi- 
photon absorption. Our theoretical analysis demonstrates that the proposed multi-subband scheme 
could represent a promising alternative to conventional quantum-well infrared photodetectors in 
the terahertz spectral region. 
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I. INTRODUCTION 



The recent development of reliable far-infrared (far-IR) semiconductor-based laser 
sources, such as the quantum-cascade (QC) laser 1,2 , together with the potential applica- 
tions in imaging, communication and medicine, identify terahertz (THz) radiation detection 
as a crucial technological milestone. To this end, many approaches have been proposed in 
the last years, which aim at accessing the 1-10 THz region of the electromagnetic spectrum. 
Currently proposed solutions encompass a variety of different approaches, each with its own 
peculiar characteristics. 

From the electronics world, field effect transistors are extending their operation fre- 
quency into the sub-THz and THz region exploiting plasmon resonance effects^^. On the 
other hand, optoelectronic techniques benefitting from electro-optical properties of LiTaOs, 
LiNbOs and ZnTe crystals have been proposed^. 

Semiconductor heterostructures also play a fundamental role in this field, QC structures^ 
as well as quantum-well infrared photodetectors (QWIPs)^ being among the most promising 
directions. Concerning the latter, radiation detection via conventional QWIP designs resorts 
on direct bound-to-continuum electronic transitions, which allowed to achieve remarkable 
levels of performance in the mid-IR range. Recently, the use of multi-level architectures, 
opening up to bound-to-bound electronic transitions, has been proposed and studied, focus- 
ing both on their intrinsic nonlinear character and on their wide-band absorption spectra. 
While the latter feature allows for multi-color— or wideband detection^^ 1 ^^, second-order 
nonlinearities of two-level systems have been studied and experimentally demonstrated with 
the idea of using the devices for second-order autocorrelation measurements.— iiLi§ii2 

The extension of the conventional, bound-to-continuum, QWIP principle into the far-IR 
range is not straightforward. In particular, one of the main issues in THz-operating devices 
are the huge dark current values that cause the background limited infrared photodetection 
temperature (Ikip) to be in the range 10-15 K,— >2I that is, much lower than that of state-of- 
the-art mid-IR QWIPs. In a previous work,— we addressed the advantages of the application 
of multi-level architecture in THz QWIP designs, and concluded that a bound-to-bound- 
to-continuum scheme may efficiently face the above-mentioned dark current issue. More 
recently^, we have analyzed the performances of such novel architecture, focusing on the 
characteristic figure of merit T blip . Our results suggest the possibility to achieve a consistent 



2 



improvement of the operation temperature of THz QWIPs by means of our proposed multi- 
level design. In the present article, our findings are further discussed and the theoretical 
model on which our calculations are based is explained in more detail. 

II. PHYSICAL SYSTEMS AND MODELING STRATEGY 

Our prototypical device consists of an infinitely periodic semiconductor-based het- 
erostructure supporting, within each period, a set of equally-spaced bound states. The 
physical system we are considering is therefore an electron gas within a periodic nanostruc- 
ture and in the presence of external electromagnetic fields. The corresponding Hamiltonian 
can be schematically written as 

H = H° + H'. (1) 

The first term of Eq. ([I]), 

H° = HI + H° qp = t«cic a + ^Jx* (2) 

a Aq 

is the sum of the free-carrier (H°) and free-quasiparticle (H° ) Hamiltonians, where the 
fermionic operator c' a (c a ) denotes creation (destruction) of a carrier in the single particle 
state a, with energy e a , while the bosonic operator (b\q) denotes creation (destruction) 
of a quasiparticle excitation of type A (phonons, photons, plasmons, etc) with wave vector 

q 

The Hamiltonian H' in Eq. (TjQ) is the sum of all possible interaction terms between 
electrons and quasiparticles. Since the aim of the present paper is to provide a focus on 
the electron-photon interaction dynamics, the latter will be treated in a fully microscopic 
scheme, in terms of the Fermi's golden rule. Conversely, all the other carrier-quasiparticle 
interactions will be described within a phenomenological electronic mean-lifetime picture, 
providing effective scattering probabilites that guarantee the proper thermalization of the 
electron population in the absence of external electromagnetic fields. 

A. Band structure calculation 

The single-particle Hamiltonian H° describes the non-interacting carrier system within 
the effective three-dimensional potential profile of our quantum device. The generic label a 
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adopted in Eq. (J2]) denotes, in general, a suitable set of discrete and/or continuous quantum 
numbers; for the case of quasi-two-dimensional semiconductor heterostructures, as the ones 
considered in this paper, the latter includes a partially discrete index along the so-called 
growth direction. In particular, since our prototypical design is made up of a sequence of 
identical units, the potential term consists, in the envelope-function formalism, of a periodic 
one-dimensional (ID) profile. 

For a device grown along the ^-direction and homogeneous as far as the in-plane (x, y) 
dynamics is concerned, the following factorization of the electron wavefunction may then be 
assumed 

^bk z k p (r) = A,k z O)0k p 0, v) , (3) 

where k p and k z are the in-plane and along- z components of the electron wavevector k, 
respectively, and b is the label numbering the various discrete subbands in which the con- 
duction band is split because of the ID quantum confinement potential. 

While parabolic bands are considered for the in-plane dispersion, and 0k p y) is the 
corresponding plane wave, the band structure along the growth direction is computed from 
the ID Schrodinger equation for the given potential profile. Moreover, due to the typically 
low doping levels in this kind of devices, charge-density effects on the potential profile may 
safely be neglected and no Schrodinger-Poisson coupling is included in our modeling. 

The Schrodinger equation projected along the z direction is solved by means of a plane- 



wave expansion, as described in Ref. 



241 ] . The following basis functions may then be adopted 



Xn ,kM = -4= e 't G »+*-> (4) 



where n is an integer running from — N to N, L z is the period of the ID potential (i.e., the 
supercell width), and G n = 2irn/L z and k z (—tt/L z < k z < tt/L z ) are the reciprocal lattice 
vector and the quasimomentum in the first Brillouin zone, respectively. The basis functions 
are normalized, as usual, over the supercell 

-L z /2 

dz = 5 nm . (5) 

L z /2 

In a reduced-zone scheme we can express the along- ,2 wave function of an electron in 
subband b and momentum k z as 

[l~ z n 

^b,kA Z ) = C b,n,k z Xn,kA Z ) ■ (6) 

n=-N 



The series expansion in Eq. ([6]) allows us to convert the stationary Schrodinger equation into 
a discrete eigenvalue problem. The solution of such a problem consists of a set of 2N + 1 
energy eigenvalues, £b,k z i each representing the allowed energy level for an electron in subband 
b with wavevector h z . The components Cb tTl ,k z represent the spectrum of the wavefunction in 
the plane wave basis set. The plane-wave like normalization of the wavefunctions ^b,k z ,k p , 

(^bk^bk^) = J Vtk> K (r)*bk^ P ( r ) dr = 5(k- k') , (7) 

is guaranteed by the form (Tj0) of the basis functions and is consistent with the fact that the 
structure is assumed to be infinite along z. 

B. Potential profile 

The quantum design of our semiconductor device should satisfy several requirements. 
First of all, the main constrain is to have equally spaced bound levels. Secondly, we want 
to be able to control the number of such levels and their spacing too. 

When speaking of equally spaced levels, the first solution would seem to be that of a 
parabolic potential profile. The implementation of the latter, however, besides non-trivial 
growth issues, poses more fundamental problems: in order to have carrier transport we need 
a continuum and thus the parabolic potential must be truncated at some point. Such a 
truncated parabola would not support equally spaced levels anymore. We therefore decide 
to use multi-quantum-well strategies for our QWIP basic period. 

Single quantum wells are used to produce the single bound level providing the bound-to- 
continuum transition exploited in conventional QWIPs. Two energetically equal transitions 
can still be obtained with a single QW of proper geometry. The tuning of the separations 
of three bound levels cannot be achieved with a potential having only two free parameters 
(width and depth) and thus we have to switch to more complex structures. 

The nested QW structures — shown in Figure turn out to be convenient choices. 
The introduction of additional geometrical parameters to the standard QW design allow 
us to control the number and position of the desired number of energy levels. A detailed 
description of the method used to determine the potential profiles is given in Appendix [A] 

Fig. [1] shows the supercells of our prototypical structures which are to be infinitely repli- 
cated along the growth direction. The use of many repetitions of the basic unit is indeed the 
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z (nm) 

FIG. 1: Potential profiles along the growth direction of our prototypical devices, designed to 
operate at 3 THz, with a number of bound states varying from one (a) to four (d). The proposed 
symmetric nested-quantum- well structure, in (c) and (d), provides additional geometric parameters, 
with respect to the single QW design, that can be varied to tune the energy level separation. 

strategy exploited in this kind of unipolar devices to optimize detection efficiency. Finite- 
size (i.e., boundary and contact) effects are therefore of minor importance. Moreover, due 
to the low doping values, in-plane quantum confinement effects are negligible too. 



C. Transport model 

The transport model we employ, to describe the electron dynamics in our unipolar device, 
is based on the Boltzmann transport equation describing the distribution of electrons in the 



6 



device conduction band. Its general form for the case of N subbands is the following 



df b (k) e 



N 

F ■ V/ 6 (k) + / [ p wfc k')^(k') - P b , b (k', k)/ 6 (k)] dk' 
b'=i J 



(8) 



dt h 

where /&(k) is the single-particle distribution function of electrons in a state with wavevector 
k in subband b, p,/&(k', k)dk' is the probability per unit time that a scattering event bringing 
an electron from a state in band b and wavevector k to a state in band b' and wavevector k' 
occurs, and F is the external electric field providing the electron drift. F may in general be 
oriented in any direction; in this paper, we will limit our discussion to biases applied only 
along the growth axis. 

The knowledge of /&(k) allows us to evaluate the current density across the device, J7, as 
follows 

^) = (^E/ V ^( k )^ k ) dk ' ( 9 ) 

provided that the distribution function is normalized as 

^^|/ 6 (k)dk = iV e , (10) 

where Pfe(k) is the miniband dispersion and N e the number of electrons per unit volume in 
the device. 

Being interested in the steady-state behavior of our device, we solve the homogeneous 
equation obtained from Eq. (jSj) when the time derivative is set equal to zero 

N . 

|F • V/ 6 (k) +J2 i p ^'(K kU'(k') - P n (k', k)/ 6 (k)] dk' = . (11) 
b>=i •' 

The latter equation is solved employing a finite difference strategy as described in Ap- 
pendix [Bl 

The various scattering mechanisms affecting the electron dynamics are included into the 
global probabilities P&/&(k',k) and may be separated into the following contributions 

P^(k', k) = P°f(k', k) + P^(k', k) , (12) 

where P^'(k', k) is the electron-photon interaction part and P 6 *^(k', k) accounts for all ther- 
malization processes. 
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D. Non-optical scattering model 



To keep the model as simple as possible, yet without spoiling the proper description of 
the main physical issues, all non optical scattering processes are accounted for by means of 
a phenomenological mean lifetime r which acts as a global fitting parameter. 

Let us introduce a thermal transition probability density P^(k, k') such that the mean 
lifetime r^k of an electron in band b with wavevector k is given by 

— = £ / ^( k '> k ) dk '- (13) 

^ J 

The mean lifetime of the electrons, r, is then defined in terms of the distribution function 
/ 6 (k) as 



r^/f*^?// fMK k) dk ' dk ' 



(14) 



66' 

tth 



The latter can be used to compute the P probabilities once r has been fixed and a functional 
form for P th has been set. However, the definition of r given in Eq. (fTl!) implies the 
knowledge of the single-particle distribution function, which is obtained from Eq. pip , which 
in turn requires P th to be determined. To break this loop we choose to drop the strict physical 
interpretation of r as the actual mean lifetime of electrons and simply use it as a measure 
of the strength of thermalization mechanisms. In this picture we can perform the mean in 
( TT4j) using a distribution function of our choice and convenience, bearing in mind that this 
won't affect our conclusions. We thus define r as 

To evaluate the integral in Eq. ( fl5|) . we have to choose a functional form for P th containing 
a free parameter suitable for normalization. Since P th must account for all thermalization 
mechanisms, its form must ensure that in the absence of any external excitation (i.e., no 
bias, no light) the system exhibits a thermal distribution function, that is, a distribution 
function such that 



AM 

At equilibrium we know from the detailed-balance principle that 

Mk') p$(k',k)' 



(17) 



the simpler way to fulfill this requirement is to impose 



yth _ 
b'k',6k — 



P V Vb (k'k) 



(18) 



Po 



1 



if E b {k) >E b ,(k') 



e 



E b ,(k')-E b (k) 



if E b (k) < E v {k!) 



where Po is a normalization constant that can be computed in terms of r as follows 



The strategy is therefore to first assume a value for r and then use the latter to compute 
Po. This completely determines the probabilities P th that appear in Eq. ( |T2l) and allows us 
to solve Eq. ffTTj) . 

Actually the definition of Po would not be of any importance if thermal scattering were 
the only scattering process, but, since we want to investigate its competition/interplay with 
carrier-photon interaction, Po (and consequently r) is the parameter that allows us to adjust 
the relative strength of the two mechanisms. 

E. Electron-photon interaction 

To evaluate the Tbn p of our prototypical device, we have to properly describe the in- 
teraction between the electron population and the radiation field of an external blackbody 
source. 

The second-quantization electric- and magnetic-field operators for a plane electromagnetic 
wave with wavevector q have the form 




(19) 




(20) 




or alternatively 




(22) 



(23) 
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where E q and B q are the classical electric and magnetic fields, u; q is the dispersion relation of 
the medium, e is the dielectric constant, \i the magnetic permittivity, V is the device volume, 
e q and b q are the polarization unit vectors such that e q ■ q = b q • q = e q ■ b q = 0, and a q 
and a q are destruction and creation operators, respectively, for a photon of wavevector q. 

The expressions above allow us to write the electric and magnetic field operators in the 
case of a linear superposition of plane waves as 

E = ^E q (24) 
q 

B = ^B q . (25) 
q 

With the latter definition, we can easily recover the usual expression for the second quan- 
tization hamiltonian of a population of photons in terms of the energy density operator 
U(r) 



H ph 



V q 



[ f)(r)dr = 
Jv 

£(^ q .Et + ^B q .Bt-_^dr= (2g) 
q 

In this picture, the classical energy density, U(r) = jS^^E^ + 2^S q -^q> re f ers to the 
zero-point energy density, , of the electromagnetic field in a cavity of volume V. 
Given the electric field operator E, we can define the vector potential operator A as 

A = E A * = E 4( A <^ + A *A) = E|f> ( 27 ) 

q q v q q 

where A q is the classical vector potential, having implicitly assumed a gauge where E = 4^. 
In a second-quantization picture, the electron-photon interaction hamiltonian operator 

Hopt = _^-V • A = -- ^= Y (V • A q a q + V • A*a q ) (28) 

m m\/2^ q q q q 

q 

can be written as 

fyopt = ^ ^ Qa , q C^a q C Q / + gaa'^a'O^Ca (29) 



aa'q 
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here, the first (second) contribution describes a process in which an electron performs a 
transition between the two single-particle states a = (6, k z , k p ) and a' = (&', k' z , k') absorbing 
(emitting) a photon; this mechanism has a coupling constant g which is expressed as 

Wq = -= / drv^(V • A q )tt Q , . (30) 

m\/2 7 

The evaluation of g QCr ' q from Eq. (1301) can be carried out in terms of the plane wave 
expansion of given in Eq.s ([3]) and (J6]) 



9bk z k p ,b'k' z K , q = -^7= 5(k + q - k') V c£, nfc , c bnfcz [A Zjq (g 2 + G n + fc z ) + A Piq ■ k p ] , (31) 
mv2 

where A z q and A p q are the along- z and in-plane components of the vector potential, re- 
spectively. 

Equation (|3T|) may be simplified in several ways. First of all, the usual dipole approxi- 
mation allows to neglect the photon momentum q with respect to the electron momentum 
k. 

Since we have assumed a parabolic in-plane dispersion, and since we expect electrons to 
have a quasi-thermal distribution, then the great majority of them will occupy states close 
to the subband bottom (k p ~ 0). On the other hand, the minibands along k z are either flat 
or slightly dispersive, that is, much narrower that the related subbands. Therefore we may 
assume that, for the majority of the electrons, k p {G n + k z ). 

This leads to the following simplified expression for the coupling constant 

SW,q= i=$(k a - k Q /)|A q | cos<y9 q )cl, nkl c bnkz (G n + k z ) (32) 

my 2 

v n 

where <y9 q is the angle between the vector potential and the z direction. The relevant term 
in the computation of transition probabilities is |5w, q | 2 which contains a cos 2 <£> q term. If 
we consider a blackbody radiation we can assume it as composed of a superposition of plane 
waves with random polarization and thus we would replace cos 2 y9 q with its mean value over 
(0,27r), that is 1/2. Anyway each electromagnetic mode is the sum of two independent 
polarizations thus we may simply replace cos 2 <^q ~ 1, obtaining 

fie 

9aa'M = 7^( k a ~ k a /)|A q |jW (33) 

my 2 
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where 

Paa' = ^ C b'nk' z C bnk z ( G n + k z ) (34) 
n 

is the matrix element of the momentum operator between states a' and a. 

Let us now consider a photon absorption process, bringing the system from state | a', n q ), 
with an electron in state a' and n photons with wavevector q, to state | a, n q — 1), with 
the electron in state a and (n — 1) photons in state q. Its probability per unit time can be 
evaluated by Fermi's golden rule as 

2ir 

P S, q = - 1 I H°*\a{,nJ\ 2 6(E a - E a , - huj . (35) 

The calculation gives 

27T 

P aa',q, = ~T I W.ql^q^ (E a ~ E a > - tlUq) . (36) 

On the other hand, the probability of a photon emission process, in which the system 
performs a transition from state | a, ra q ) to state | a', n q + 1), is 

2-7T 

P Xq = y k'cqi'K + l)5(^a ~ + ^ q ) . (37) 

F. Interaction with blackbody radiation 

Since our aim is to determine the T b i ip of our prototypical detector, we need to study its 
interaction with the background radiation, considered as a blackbody radiation at 300 K. 
From a quantum mechanical point of view, a blackbody radiation is a photon population at 
thermal equilibrium following the Bose-Einstein distribution law. 

A non interacting electron system only coupled to a photon bath at thermal equilibrium, 
must itself thermalize. Indeed, by employing the detailed-balance principle and substituting 
the Bose-Einstein distribution in Eqns. (I3T)|) and (|37[) we can write 

f p°Pt fi B a -E_j 

J a aa',q _ "■q _ — v — w 



fa' P°f (riq+1) 
a'a,q 



e k B T , (38) 



that is, the steady-state distribution function is such that the ratio between the occupation 
numbers of states a and a' is, as expected, the Boltzmann factor. 

Equations (|36|) and (l37j) give the transition probabilities for an electron interacting with 
an electromagnetic plane wave, which can be seen as an electromagnetic mode of a cavity. 
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When our device is inside a cavity at thermal equilibrium (a blackbody), the total transition 
probabilities must be summed over all modes q. This is also formally described by the 
interaction hamiltonian (1291) which is a sum over all wavevectors q. We therefore write, for 
the absorption process, 

P£ = E P -', q =yE lw, q |W(#« " E a , - huj . (39) 



% 

q q 



2 



In the limit of a infinitely large cavity, the summation becomes an integral in dq and \g a 'a,m 
becomes a spectral density |^a'a(q)| 2 which is related to the squared vector potential spectral 
density |^4(q)| 2 through equation fl33|) . 

The quantity |A(q)| 2 can be expressed in terms of the energy density U(q) = |e|E(q)| 2 + 
2^|B(q)| 2 and, considering the relations E(q) = zo; q A(q), B = V x A, as 

|A(q)| 2 = ^. (40) 

£o; 2 

For a linear dispersion relation u; q = cq, we can switch from spectral densities in the wavevec- 

tor domain to spectral densities in the frequency domain. In particular, we can consider the 
spectral energy density 

,3 

- (41) 

which is the energy density of an infinite cavity in which each mode is populated by one 
photon. The term T is a constant expressing the limited field-of-view (FOV) of the device 
and depending on how the blackbody radiation is coupled into the detector in the specific 
experimental setup. 

The absorption probability can then be evaluated in the frequency domain as 

P% = fr / IwMlV")* (Aoo -u)duj= (42) 
= -^2-|w(A^)| 2 ™(ZW) 

where Au = — h a ' is the resonance frequency of the transition. 

After substitution of equations (1331) and (j4"0]) into Eq. (H2l) . we obtain 

PS = U(Au)-^\p aa ,\M^)6(k a - k a r) . (43) 

Since we are dealing with a thermal population of photons we take n(Au) as the Bose- 
Einstein distribution function so that the total absorption probability can be finally written 
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as 



ATT 2 c 3 m 2 e 

Analogously, the emission probability is 

pgpt _ e 2 hAuF 2 

G. Fixing the value of r 



P Z = ,_ 2 ^ 2 J iWr^Ai 5(k Q - k a ,) . (44) 

e k B T _ I 



e k B T _ I 



+ 1 



5(K' ~ K) ■ (45) 



The model contains a free parameter, r, which has to be adjusted in order to reproduce 
some experimental data. Its value is, in principle, crucial in determining the Tbu p of the 
simulated devices, since changing the mean lifetime of electrons will change the strength of 
thermal scattering with respect to optical scattering and thus will affect the point at which 
these two competing processes balance. 

In particular, we choose to adjust r in order to reproduce the measured T b i ip (12 K) of 



the bound-to-continuum QWIP operating at 3.2 THz and reported in Ref. [21]. 

Figure [2] shows the total normalized current densities that we obtain for a 3 THz QWIP 
as a function of temperature for different values of r. Although r is the key parameter that 
fixes the value of the Tbn p , it can be noted from the figure that in the interval r m 50 — 100 
ps the Tbnp shows little variation around 12 K. We can thus safely assume for r any value 
in this range, like, e.g., r = 80 ps, in order to reproduce the experimental data. 

It is important to stress once more, at this point, that this very large value derives from 
the fact that we are using a simplified model for thermal scattering; our fitting parameter r is 
not to be taken as a realistic indication of electron scattering time in the real heterostructure. 

Once the value of r has been set, on the basis of the above discussion, we use it in 
modeling the current response of detectors operating at identical frequencies but employing 
the proposed bound-to-bound-to-continuum strategy, and differing in the number of bound 
states. 



III. RESULTS AND DISCUSSION 

We now apply the model of Section II to describe four different devices, having a number 
of bound levels ranging from one (standard QWIP) to four, and designed according to 
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FIG. 2: Estimated normalized current density along the growth direction, as a function of device 
temperature, for the one- level QWIP of Figure 1, operating at 3 THz in the presence of a back- 
ground radiation field at 300 K. Different symbols correspond to different values of r. Applied 
electric field is 50 V/cm, FOV is 90°. The dashed line marks the current doubling. 

our bound-to-bound-to-continuum strategy. All devices are exposed to a 300 K blackbody 
radiation under a 90° FOV and are subject to a 50 V/cm external bias. 

Figure [3] shows the total normalized currents across each of the four devices as a function 
of the device temperature. Each curve allows to identify a low-temperature regime in which 
the dark current is negligible with respect to the photocurrent: the total current is therefore 
independent from the device temperature. Conversely, in the high temperature region, 
the dark current increases almost exponentially so that the photocurrent quickly becomes 
negligible and the current is totally due to the 'dark' contribution. The Tbn p may be identified 
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FIG. 3: Estimated current density (normalized) along the growth direction, as a function of device 
temperature, for the four diverse multilevel designs of Figure 1, differing in the number of bound 
states. 

as the temperature at which the total current doubles with respect to the low temperature 
region (dashed horizontal line in Fig. [3]): at this temperature the dark current and the 
photocurrent have the same magnitude. 

The diverse designs have values of Tbn p = 11.5, 19.5, 23.5 and 28.5 K for one, two, three 
and four bound levels, respectively, showing the trend reported in Figure HI 

The increase of TbH P may be better interpreted by looking at Figure [5j where the pho- 
tocurrent and the dark current are plotted as a function of the number of bound levels. Both 
currents decrease on increasing the latter, but the dark current does it faster. Therefore, 
the temperature at which the two are equal moves towards higher values. 
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FIG. 4: Estimated values of Tbn p for the four devices of Figure 1, as deduced from the data shown 
in Figure O The dashed line is a guide to the eye. 

As can be observed in Fig. [5] there is a dramatic decrease in the photocurrent when 
switching from two to three bound states, which is mainly due to the reduction of the 
photoconductive gain. In fact, for the four-level design, the latter reduces to just the 0.2% 
of the value of the one-level QWIP. Conversely the quantum efficiency is only lowered by 
14% and thus its variation does not significantly affect the photocurrent. 

This behavior can be explained by considering Fig. [Hand noting that there is a remarkable 
geometrical difference between the two-level and the three-level design. The presence of 
the nested QW introduces a new ground state whose wavefunction has little overlap with 
the wavefunctions of higher energy states, therefore reducing the oscillator strength. This 
conclusion is supported by the fact that the photocurrent reduction between the three and 
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FIG. 5: Photocurrent (squares) and dark current at 17 K (circles) as a function of the number of 
levels in the device. The dark current decreases more rapidly than the photocurrent, on increasing 
the number of bound subbands, and thus the 7bii p increases. Lines are a guide to the eye. 

four levels designs, where no major structural change has been introduced, is comparable to 
the decrease between one and two levels cases. In the present work the optimization of device 
performance is not the central issue; in this respect, a more elaborate tuning of the device 
geometry to achieve higher oscillator strengths would surely allow for better operational 
results. 
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IV. SUMMARY AND CONCLUSIONS 



The scaling-down of QWIPs to access the terahertz range of the electromagnetic spectrum 
is not straightforward: in this frequency range the dark current, mainly due to the high- 
energy tail of the electron distribution function, may indeed become predominant over the 
photocurrent signal. In a recent paper— we have proposed and theoretically investigated a 
THz-detector design alternative to the conventional QWIP structure. The former, instead 
of resorting on the conventional bound-to-continuum scheme, exploits a bound-to-bound- 
to-continuum strategy. In particular, a ladder of equally-spaced bound levels is employed, 
whose energy-step is tuned to the desired detection frequency. 

Our previous analysis demonstrated that a multilevel architecture can indeed satisfac- 
torily face the dark-current problem in far-infrared QWIPs. In the present paper we have 
significantly improved some features of our model to better reproduce the behavior of realis- 
tic state-of-the art designs. In particular, our attention has been devoted to a specific figure 
of merit of QWIPs, such as the background-limited infrared photodetection temperature 
(Tbiip), which is related to the interplay between dark current and photocurrent. Our results 
have demonstrated that the proposed multi-subband scheme allows for higher T b i ip values, 
with respect to conventional QWIP designs operating at the same frequency, and therefore 
could represent a better alternative for THz radiation detection. 



APPENDIX A: POTENTIAL PROFILE CALCULATION 

The key-point of the proposed bound-to-bound-to-continuum architecture is the design 
of the nanostructure potential profile. The latter requires the solution of the inverse problem 
of setting the desired energy spectrum and then finding the corresponding operator, i.e., the 
potential energy term of the electron hamiltonian H°. In the present paper, this problem 
has been solved numerically, by means of a variational approach. 

In particular, starting from the function V(z) that describes the potential profile in one 
period of our device, a functional ^ r [l / (z)] has been defined whose value represents how 
far the function V(z) is from our target function V(z). The latter must be such that the 
operator 

h 2 d 2 
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has a spectrum composed of a lower part with N equally-spaced discrete values (bound 
states) and an upper continuous part. The number N and the energy spacing E of the 
bound states are our design constraints. 

In general, a function V(z) will produce N bound states of energies Ei, i = 1, . . ., N, 
where N can range from one to infinity (infinitely deep potential well). If there are two or 
more bound states we can define a mean interlevel spacing 

1 N 

i=i 

and a level spreading 

1 N 

a = N^ {Ei ~ Ef ■ (A3) 
i=i 

In terms of the quantities E and a specified above, the functional .F[y(z)] is defined as 

T[V{z)\ = (1 - 5 Nf} ) +a + \E-E\. (A4) 

From the latter equation we see that T is always positive and assumes the minimum, null, 
value only when N = N, a = and E — E. 

The actual existence of the minimum depends on the functional space we choose for V(z). 
Indeed, we already know a solution for the problem T = which is the harmonic oscillator, 
but the latter cannot be taken into consideration because it is not a realistic potential profile 
and its spectrum does not contain a continuous part. 

Without going into the rigorous mathematical definition of the space, which is beyond the 
scope of the present paper, V(z) must be a periodic function with period L z . We therefore 
define V(z) on the domain —L z /2 < z < L z /2 and impose V(—L z /2) = V(L z /2). Of course 
we don't want V(z) to diverge at any point and in addition we want it to be as close as 
possible to realistic and technologically accessible potential profiles. 

We choose to take V(z) piecewise constant on its domain so that it can be described by 
a discrete set of M parameters representing widths and depths of every constant sector. In 
this way, V(z) can be represented by a point in an M-dimensional space and T actually 
becomes a function of M variables. In practice, V(z) takes the form of a multi quantum 
well or a nested quantum well structure, in which we vary depths and widths of the diverse 
layers. 
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The minimization of T is not trivial mainly because of the S N fj term that makes it 
discontinuous in an unpredictable way: by slowly varying the free parameters, the potential 
profile can suddenly produce a new bound state or lose one causing T to jump by ±1 and 
thus we cannot use methods that seek a local minimum following the function gradient. 

We adopted the easiest possible solution: starting from an initial guess for V(z), the free 
parameters are varied within a certain range to see if a minimum is present and whether 
the latter is actually the absolute one, for which JF = 0. The existence of such minima 
mostly depends on the number M of free parameters that can be varied. We choose to start 
with the minimum number of parameters (which is two for a single quantum well) and then 
gradually increase this number in order to generate more bound states. 

APPENDIX B: BOLTZMANN EQUATION SOLUTION 
1. State space discretization 

The electron dynamics in our prototypical quantum device is described by the Boltzmann 
transport equation (fTTT) . The latter will be solved by finite difference discretization of the 
derivatives and Reimann discretization of the integral. 

Due to the cylindrical symmetry of the physical problem, guaranteed, as in our case, by an 
external field applied along the growth direction only, a convenient starting point is to employ 
cylindrical coordinates, with k z being the perpendicular (growth direction) wavevector, and 
k p and 9 the modulus and anomaly, respectively, of the in-plane wavevector. 

A central difference approximation of the derivatives along k z may then be applied, with 
periodic boundary conditions accounting for the repetition of the Brillouin zone. In partic- 
ular, the values of k z span the first Brillouin zone — ^-<k z <^- forming a uniform grid of 
step Nu. The width of each discrete cell is therefore equal to Ak z = , 2 Z . 

The in-plane angle 9 is uniformly discretized in the domain [0, 2ir). The number of discrete 
cells is Ng and their size is A9 = j^. Again, central difference approximation of derivatives 
and periodic boundary conditions are adopted. 

Discretization along k p poses the problem of limiting the in-plane k-space. The electron 
distribution function f(k p ) at thermal equilibrium has the form 




(Bl) 
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and thus decays rather quickly in k p . We expect the electron non-equilibrium distribution 
to decay more or less in the same way in the presence of the external radiation field. We 
therefore set a cutoff value, / cut , below which f(k p ) is considered to be negligible, and use 
it to define a maximum value for k p in the following way 

*r = /^V^ • (B2) 

The interval [0, &™ ax ] is then discretized into a uniform grid of dimension N^ pl and derivatives 
are approximated by central difference formulae in the inner nodes. In particular, f(k p = 0) 
is supposed to have null derivative (gaussian-like behavior) and the same applies for /(A;™ ax ). 

Lmax 

Discrete cells along k p have a width Ak„ = -ft — . 

To complete the description of the state-space we need to set the number N b of subbands 
actually considered for calculations. The plane-wave solution of the Schrodinger equation 
requires from one to two hundred plane waves to give stable energy level values, producing the 
same number of subbands. However, the electron distribution function decays rather quickly 
and a number of bands from five to a few tens is usually enough to ensure convergence. 
Indeed, the actual number of bands depends on the operating conditions of the device, such 
as the temperature and the presence of incident light or external bias. 



2. Discrete Boltzmann equation 

After discretization of the state-space, the distribution function /fi(k) can be itself dis- 
cretized into a vector of components /j. The label i ranges from 1 to the total number of 
grid points N = N^N^NgN^ and accounts for the band index \ and the three k-space 
coordinate indexes k Zt i, k P) i and 9^, collectively named k;. The value /, is the mean of /^ (kj) 
over the grid volume element Ak, ; = Ak z ^Ak Py iA9i 

fi = ^~l A,(ki)dk (B3) 

from fi we define the occupation number rii of the i— th discrete cell as 

ni = /iAki . (B4) 

With this definition and taking into account relation (imp the occupation number is normal- 
ized as 

J2 n * = N e- (B5) 

i 
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After the discretization of the distribution function, we need to find a suitable discretiza- 
tion of the scattering probabilities P\,y (k, k') . The total number of particles Rji that perform 
a transition from the volume Akj in band 6j to the volume Akj in band bj, is given by the 
probability that a particle in kj performs a transition towards one of the states in volume 
Akj, which is J Ak P b:jbt (k, kj)dk, integrated over all states of the starting volume 

Rji= ( I AXk')A A (k,k')dkdk'. (B6) 

In the discretized system the probability that a particle performs the same transition is 
PjiAkj and the number of particles in the starting volume is rij = /jAkj, thus 

Rji = fiAkiPjiAk, . (B7) 



The combination of equations flB6j) and flB7j) allow us to derive the following expression 



for Pji 



If we approximate /^(k') with its mean value fi over the volume Akj, we can take it out of 
the integral and rewrite Eq. flB8j) as 

Wji = PjiAkj = 4r I I (k, k')dkdk' (B9) 

Ak i JAki J Akj 

and Eq. (JH2|) as 

Rji = WjiUi . (BIO) 

The quantity ~Wj% is the probability that an event bringing an electron from a state % to one 
of the states in volume Akj occurs. By multiplying Wji by the number of particles in the 
volume Akj, that is rii, one obtains the rate Rji of particles leaving the volume Akj and 
entering the volume Akj. 

For the generic z-th volume element we can then write a rate equation in the usual form 

= YliWijnj - WjiUi) (Bll) 

3 

which is the Boltzmann equation for a discrete system composed of N states. For the 
stationary state we write 

Y^iWijUj - Wjim) = . (B12) 
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Due to the discretization procedure, the drift term in Eq. (jSJ) can be written as 

V/ 6 (k)-|F = J]Wf/, (B13) 

3 

where is an equivalent scattering matrix that can be included into the Wij term in 
Eq. (1B12I) . The latter will in general consist of several contributions 

w a = W S + w§ + Wf (B14) 

where W$ and W%* are the discretized probability densities P^ a , and P°^' defined in 
Eqn.s ( TTHI) and ( l44l) and computed using Eq. ( 1B9I) . 



3. Non-optical scattering probabilities 

The derivation of Section Hi Dl can be followed in the discretized system by replacing /t(k) 
with fi and P&&'(k, k') with P^-. In particular, the discretized version of Eq. (JTSJ) 

^=(^E^Ak ( (B15) 
leads to a definition of P/j 1 similar to Eq. ffTS]) 

I 1 if P 7 > ^ 

P^ = P V id = P \ ^ (B16) 

I e fc s T if Pj < Pi 
where again P is a normalization constant that can be computed as 

T = WfN^ Vi >^- (B17) 

The discrete thermal transition probabilities are then written according to Eq ( 1B9I) 

w th = p th Aki _ (Blg) 



4. Optical scattering probabilities 

The discrete optical scattering probabilities can be directly computed using Eqns. flB9j) 
and (j44p or (j45j) . In both cases, the transition probabilities are of the form 

Pf (k i ,k J )=w lJ S(k l -k 1 ) (B19) 
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where Wij contains all the coefficients and differs for absorption and emission processes. 
Substituting the latter expression into Eq. flB9j) gives 

Wf = ^-\ [ Wij 5(k - k')dk'dk (B20) 



*■? Ale . 

JAkj J Alt, 



Wijdk . 



Assuming that Wij is a smooth function over the volume cell Ak,- and for a sufficiently 
dense grid, we can approximate constant and take it out of the integral, which in 

turn results to be Ak,-, therefore yielding 

Wf = , (B21) 

which is the discrete transition probability to be used in the solution of the discrete Boltz- 
mann equation. 
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